Instructions

Select between five and 30 locations within a city or metropolitan area. For each location, generate isochrones for the same travel time for at least two different modes. Calculate the area of each isochrone and compare the areas of the isochrones for the two (or more) modes you analyzed. Create three figures (which may or may not be maps) to illustrate the results of your analysis.

Setup

library(osmdata)
library(opentripplanner)
library(tidyverse)
library(sf)
library(ggthemes)
library(ggspatial)
library(ggmap)
library(stringr)
library(sp)
library(rgeos)
library(tidygeocoder)
library(raster)
childcareaddress <- read.csv("2miradiuschildcare.csv")
childcareaddress <- as.data.frame(childcareaddress)
opq(bbox = 'Boston MA USA') %>%
  add_osm_feature(key = 'highway') %>%
  osmdata_xml(file = 'OTP/graphs/default/boston_streets.osm')
MA_state_plane <- "+proj=lcc +lat_1=41.71666666666667 +lat_2=42.68333333333333 +lat_0=41 +lon_0=-71.5 +x_0=200000 +y_0=750000 +ellps=GRS80 +units=m +no_defs"

boston_street_features <- opq(bbox = 'Boston MA USA') %>%
  add_osm_feature(key = 'highway') %>%
  osmdata_sf()

boston_streets <- boston_street_features$osm_lines %>%
  st_transform(crs = MA_state_plane)
path_otp <- otp_dl_jar("OTP")
otp_setup(otp = path_otp, dir = path_data, memory =1024)
otpcon <- otp_connect()
## Router http://localhost:8080/otp/routers/default exists

Geocoding Addresses of Childcare Centers in Mattapan

I found this list of addresses of certified childcare centers in Mattapan from the Mass.gov website. I cleaned up the data in excel to comport with geocoding the addresses using tidygeocoder.

address_list = c("35 Astoria Street, Mattapan, MA",
"40 HOSMER STREET, Mattapan, MA",
"32 WOODBOLE AVENUE, Mattapan, MA",
"7 LANDOR ROAD # 1, Mattapan, MA",
"23 WILMORE STREET #1, Mattapan, MA",
"5 Mildred Ave, Mattapan, MA",
"500 Walk Hill Street, Mattapan, MA",
"97 WEST SELDEN ST, Mattapan, MA",
"608 NORFOLK STREET, Mattapan, MA",
"130 Orlando Street, Mattapan, MA",
"62 RIDGEVIEW AVE, Mattapan, MA",
"55 Hollowell Street, Mattapan, MA",
"760 Morton Street #5, Mattapan, MA",
"100 ROCKDALE ST, Mattapan, MA",
"10  ELLISON AVE, Mattapan, MA",
"6 Lorna Rd, Mattapan, MA",
"1295 BLUE HILL AVE, Mattapan, MA",
"14 Groveland Street #1, Mattapan, MA",
"70 FAVRE ST # 2, Mattapan, MA",
"124 LORNA ROAD, Mattapan, MA",
"37 Duke Street  #2, Mattapan, MA",
"14 DANIA STREET, Mattapan, MA",
"37 Babson Street, Mattapan, MA",
"5 Mildred Ave, Mattapan, MA",
"778 MORTON STREET # 1, Mattapan, MA",
"624 HARVARD STREET, Mattapan, MA",
"438 River Street, Mattapan, MA",
"535 RIVER STREET, Mattapan, MA",
"8 LESTON ST, Mattapan, MA",
"596 RIVER ST, Mattapan, MA",
"5 MILDRED AVENUE, Mattapan, MA",
"22 HARMON STREET, Mattapan, MA",
"73 WOODBOLE AVENUE #216, Mattapan, MA",
"74 WESTMORE ROAD, Mattapan, MA",
"5 Mildred Avenue, Mattapan, MA",
"712 MORTON STREET #3, Mattapan, MA",
"119 Lorna Road, Mattapan, MA",
"36 Rockingham Road, Mattapan, MA",
"29 Woodruff Way, Mattapan, MA",
"9 GREENDALE RD # 2, Mattapan, MA",
"54 Hiawataha Road, Mattapan, MA",
"1009 Morton Street, Mattapan, MA",
"163 COLORADO ST, Mattapan, MA",
"255 RIVER ST, Mattapan, MA",
"152 STANDARD STREET # 129, Mattapan, MA",
"55 Woodhaven Street Apt. 2, Mattapan, MA",
"12 ELLISON AVENUE # 1, Mattapan, MA",
"557 NORFOLK ST, Mattapan, MA",
"100 Hebron Street, Mattapan, MA",
"20 Outlook Road, Mattapan, MA")
points <- geo(address = address_list, mode = "batch") 
head(points)
## # A tibble: 6 x 3
##   address                              lat  long
##   <chr>                              <dbl> <dbl>
## 1 35 Astoria Street, Mattapan, MA     42.3 -71.1
## 2 40 HOSMER STREET, Mattapan, MA      42.3 -71.1
## 3 32 WOODBOLE AVENUE, Mattapan, MA    42.3 -71.1
## 4 7 LANDOR ROAD # 1, Mattapan, MA     42.3 -71.1
## 5 23 WILMORE STREET #1, Mattapan, MA  42.3 -71.1
## 6 5 Mildred Ave, Mattapan, MA         42.3 -71.1
points2 <- na.omit(points)
points2 <- st_as_sf(x = points2,                         
           coords = c("long", "lat"),
           crs = 4326)

I then re-merged the geocoded addreses with the original csv in order to not lose the data associated with the addresses.

merged.data <- merge(points2, childcareaddress, by="address")

Calculating Isochrones

With my geocoded childcare addresses ready to go, I calculated the 5 minute walking and biking isochrones around each address.

multiple_points_5min_walk <- 
  otp_isochrone(otpcon = otpcon, fromPlace = points2, 
                mode = "WALK", cutoffSec = 300)
multiple_points_5min_bike <- 
  otp_isochrone(otpcon = otpcon, fromPlace = points2, 
                mode = "BICYCLE", cutoffSec = 300)

Map of 5-minute Walkshed around Mattapan Childcare Centers

ggplot(multiple_points_5min_walk) +
  annotation_map_tile(zoomin = 1, progress = "none") +
  geom_sf(fill ="blue", alpha=0.2) +
  theme_map() 
## although coordinates are longitude/latitude, st_intersects assumes that they are planar
## although coordinates are longitude/latitude, st_intersects assumes that they are planar

Map of Mattapan Childcare Center 5-Minute Walksheds by Access to Bus Stops

bus_stops <- st_read("https://opendata.arcgis.com/datasets/55586c8f54954f8e8fae5f40cb953d15_0.kml?outSR=%7B%22latestWkid%22%3A26986%2C%22wkid%22%3A26986%7D",
                     quiet=TRUE)
multiple_points_5min_walk <- multiple_points_5min_walk %>%
  mutate(transit_score = lengths(st_covers(geometry, bus_stops)))
ggplot(multiple_points_5min_walk) +
  annotation_map_tile(zoomin = 1, progress = "none") +
  geom_sf(aes(fill=transit_score), alpha=.5) +
  theme_map() 
## although coordinates are longitude/latitude, st_intersects assumes that they are planar
## although coordinates are longitude/latitude, st_intersects assumes that they are planar

Map of Mattapan Childcare Center Isochrones by 5-min Walk- and Bike-shed

iso_5min_walk <- 
  otp_isochrone(otpcon = otpcon, fromPlace = points2, 
                mode = "WALK", cutoffSec = 300) %>%
  st_transform(crs = MA_state_plane) %>%
  mutate(mode = "walk")

iso_5min_bike <- 
  otp_isochrone(otpcon = otpcon, fromPlace = points2, 
                mode = "BICYCLE", cutoffSec = 300) %>%
  st_transform(crs = MA_state_plane) %>%
  mutate(mode = "bicycle")

iso_all_modes <- rbind(iso_5min_walk, iso_5min_bike)
ggplot(iso_all_modes) +
  annotation_map_tile(zoomin = 2, type = "cartolight", progress = "none") +
  geom_sf(aes(fill = mode), alpha = 0.5) +
  geom_sf(data = points2) +
  scale_fill_viridis_d(name = "Area that is reachable within 5 minutes",
                       labels = c("By bike", "By foot")) +
  theme_map() +
  labs(caption = "Basemap Copyright OpenStreetMap contributors")

Scatterplot of Mattapan Childcare Center Isochrones by 5-min Walk- and Bike-shed

iso_areas <- iso_all_modes %>%
  mutate(area = st_area(iso_all_modes)) %>%
  st_set_geometry(NULL) %>%
  pivot_wider(names_from = mode, values_from = area) 
iso_areas <- iso_areas %>%
  filter(bicycle != "NULL", walk != "NULL") %>%
  filter(str_detect(bicycle,"c")==FALSE)
ggplot(iso_areas, 
       aes(x = as.numeric(walk), y = as.numeric(bicycle))) +
  geom_point() +
  scale_x_continuous(name = 
            "Area within a five-minute walking distance\nof a childcare center\n(square km)",
            breaks = breaks <- seq(10000, 260000, by = 20000),
            labels = breaks / 1000000) +
  scale_y_continuous(name = 
            "Area within a five-minute driving distance\nof a childcare center\n(square km)",
            breaks = breaks <- seq(0, 2000000, by = 100000),
            labels = breaks / 1000000) +
  theme_bw()

Scatterplot of Mattapan Childcare Centers by 5-minute Walk- and Biked-Shed and Center Capacity (Number of Children)

mergedgeometry <- read.csv("mergeddatafromPlace.csv")
mergecapacityiso <- merge(mergedgeometry, iso_areas, by="fromPlace")
ggplot(mergecapacityiso, 
       aes(x = as.numeric(walk), y = as.numeric(bicycle), size = Capacity)) +
  geom_point() +
  scale_x_continuous(name = 
            "Area within a five-minute walking distance\nof a childcare center\n(square km)",
            breaks = breaks <- seq(10000, 260000, by = 20000),
            labels = breaks / 1000000) +
  scale_y_continuous(name = 
            "Area within a five-minute biking distance\nof a childcare center\n(square km)",
            breaks = breaks <- seq(0, 2000000, by = 100000),
            labels = breaks / 1000000) +
  scale_size_continuous(name = "Childcare Center Capacity")+
  theme_bw()